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1. Introduction 

Initial State Radiation (ISR) is the most important QED correction to the Born cross 
section. For instance, in the case of electron positron annihilations at energies beyond 
the Z-threshold it leads to the so-called radiative return causing a hugely increased cross 
section. There, via multiple emission of photons, the electron positron pair tends to reside 
on the Z-resonance (returning on it) picking up the propagator in its resonance region. 
In this way, the correction term becomes several times as large as the Born term at the 
initial center of mass energy. Evidently, the effect of ISR is of great significance for the 
interpretation of experimental results. However, since initial state photons tend to be 
collinear and soft they easily escape detection. In the detector, they tend to leave only 
imprints by a reduced overall energy and some netto transverse momentum and - possibly 
- by a boost of the visible final state along the beam axis. Due to this rather indirect 
method of measuring it, the precise simulation of ISR is crucial for most experimental 
analyses. 

Several strategies exist for this task: 

1. The structure function approach 

The idea behind this approach is to simulate ISR by a probability density to find 
an electron with reduced momentum inside an incoming electron. Accordingly, the 
missing momentum is then attributed to photon radiation. That way the total cross 
section <r tot including ISR is defined as the leading order cross section ctb convoluted 
with two structure functions T ee : 

0tot = dxidx 2 T ee (xi,s)T ee (x2,s)aB{xiX2s) . (1.1) 
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Within the structure functions logarithmically enhanced contributions are exponenti- 
ated and thus re-summed to all orders with additional terms covering non-factorisable 
QED effects up to a finite order a n . Additionally, terms describing weak and strong 
effects can be taken into account. However, since the structure function is integrated 
over all photon energies and transverse momenta, it is limited to situations in which 
only the overall effects of energy reduction and longitudinal boost via ISR are of 
interest. Arbitrary phase space cuts are difficult to implement. Specific photonic 
observables like for instance the number of photons above some energy threshold or 
their transverse momenta are beyond the reach of this approach. 

2. The parton shower approach 

In this approach, the electron radiation is followed step by step from the comparably 
low scales of the order of the electron mass connected to the incoming electrons 
to the high energy scale of the actual hard process. After setting up an electron 
distribution function at the initial scale, any individual photon emission is governed 
by the well-known splitting function 

leading to a Dokshitzer-Gribov-Lipatov-Altarelli-Parisi (DGLAP) evolution equa- 
tion for the structure function. Technically, the original structure function is being 
reproduced statistically at each scale via backward evolution and re-weighting, where 
the analytic integration is replaced by Monte Carlo integration. Apparently, photons 
are generated explicitly which allows for instance the application of arbitrary phase 
space cuts. However, one of the weaknesses of this approach is that it re-sums only 
the logarithmic enhanced contributions in the DGLAP equation. Even though the 
radiation of the hardest photon can be corrected with the first order matrix element, 
it is not clear how to systematically improve to higher perturbative order, i.e. how 
to make sure that matrix element expressions for the emission of additional photons 
are recovered for arbitrary processes and arbitrary numbers of photons. 

3. The approach by Yennie, Frautschi and Suura (YFS) [1] 

Within this approach the explicit generation of photons can be corrected systemati- 
cally to all orders in the coupling constant. The basic idea here is to introduce some 
arbitrary infrared cut-off on the photon energies and treat low energy real photons 
as "un-resolvable" . Their respective contribution will cancel the emerging virtual 
infrared divergences to all orders leading to finite factorisable terms which can be 
exponentiated. In this framework correction terms for harder photons can be easily 
introduced. Since this is the approach we are going to follow we postpone a more 
thorough discussion. 

In this paper we want to give a description of the YFS method implemented in AMEGIC++ 
[2]. The outline is as follows: In section 2, after briefly reviewing the YFS approach and 
a short introduction to the nomenclature, we will introduce our new approach to calculate 
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matrix element corrections, where we make use of the ability of AMEGIC++ to calculate 
(almost) arbitrary tree level matrix elements automatically. In fact, it can be anticipated 
that our method might allow for the construction of a more generic simulation code for 
QED initial state radiation in arbitrary processes in the near future. In section 3 we 
will discuss the actual implementation of our approach and the interplay of ISR with the 
calculations of matrix elements and cross sections. We justify our method in section 4 by 
confronting our results for the radiative return to the Z-pole with experimental data from 
the LEP experiment. The generality of our program is illustrated by the application to the 
s-channel production of Higgs bosons at a possible muon collider. We will conclude with 
summarizing remarks. 



2. Method 



In this section we want to discuss our new method to include the effect of exact matrix 
elements in the generation process of initial state photons along the lines of the YFS 
approach. We will begin with a mini-review of this method suitable for Monte Carlo 
implementation. 

Mini-review of the YFS approach 

The idea underlying the YFS approach is to separate the phase space for the emission of 
real photons into two regions via a cut-off e on the energy fractions, such that photons are 
coined infrared if their energies w 7 < eEbcam • The contribution of these infrared photons is 
then used to cancel the virtual infrared divergences order by order in a. The remainders of 
this procedure factorize and can be exponentiated into an universal factor, the YFS form 
factor: 



F YFS (e) = exp 
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The functions S(k) and S(k) denote the well known universal factorizing "radiation factors" 
or "eikonal factors" for virtual and real photons for a pair of two external charged lines 
with four-momenta p± and p2 
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The emission of the visible real photons can then be corrected systematically to all orders 
in a to reproduce exact results as given by the corresponding matrix elements. To be 
more specific, let us consider as an example 2^2 processes of the type l + {pi)l~{p2) — ► 
f(qi) f ((fc), where the pi and qj label the four-momenta of the external particles. The total 
cross section including arbitrary numbers of real or virtual photons can be written as 



Y[-nr s(pi,P2,ki) e 




— £ 



(2.5) 



The YFS form factor F YFS (e) covers the contribution of factorizing soft real and virtual 
photons to all orders. The integral over the phase space of the final state particles consists of 
the integral over the two outgoing momenta q\ and q<i plus a sum over all possible numbers 
of photons ki with energy fractions above the resolution threshold e. This constraint is 
reflected in the ©-functions. The conservation of total four-momentum is enforced by the 
c)-function in the second line. Finally, the last line includes the matrix element corrections. 
The (3 n denote infrared safe combinations of cross sections and eikonal factors with n 
additional real photons as well as all finite contributions of any number of virtual photon 
loops. In practice, the order of virtual contributions is limited. In (3$ the superscript I 
denotes the total number of virtual and real photons taken into account, i.e. the order of 
perturbation theory. From there we can read off the number of virtual photons, given by 
V-n). 

For a treatment exact up to first order a, it is sufficient to determine fi^ , fi^ and f5^ , i.e. 
the cross sections for l + (pi)l~{p2) f(Qi)f(l2) up to one loop and for l + (pi) l~ (3)2) ~^ 
/fe) l(ki) at the tree level. They are defined as 

P[ 1 \k i ) = p?(k i )-S(k i )4° ) • (2.6) 

The Pr) are the fully differential cross sections for the hard process l + (pi)l (P2) ~ * 
f(qi) /(<?2) with n real photons added and up to Ith order in a, and S denotes the eikonal 
factor as defined in Eq.(2.4). 

Analyzing the master equation, Eq.(2.5), we realize that 

1. the number of photons that are explicitly produced is not constrained, 

2. the resolved photons are allowed to have transverse momentum, 

3. the cross section is independent of e (in fact, only the number of explicitly produced 
photons depends on the resolution scale), and 
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4. that in contrast to the parton shower approach any photon emission is corrected by 
the appropriate matrix element (s) to the chosen order in a including interference 
effects in multi photon final states. 

The MC algorithm 

A Monte Carlo implementation for the YFS approach outlined above was described in [3] 
for the first time. The proposed procedure there was 

1. to chose the cm. energy with the form factor F YFS , 

2. to determine the number of resolvable photons, 

3. to supply each of the photons with a four-momentum k according to the eikonal 
factor S(k), 

4. to reject all momenta outside their kinematical limits, and 

5. to correct the configuration chosen with an additional weight stemming from the 
exact matrix element. 

A more detailed description of the Monte Carlo steps can be found in Appendix A. 
Implementing matrix element corrections 

The corrections on the photon distributions and the cross sections due to exact matrix 
elements can be read off from Eq.(2.5). Up to first order in a the corresponding weight 



where again the pi are the incoming momenta, the qj are the outgoing momenta and the k n 
are the momenta of the explicitly generated photons. In the following we will again restrict 



and two outgoing momenta only. However, the extension to more final state particle is 
straightforward. 

An implementation of this weight might seem trivial, however, one serious question re- 
mains: How do we calculate matrix elements in situations, where the number of generated 
photons exceeds the number incorporated in the matrix element? Matrix elements are only 
unambiguously defined if all incorporated momenta comply with four-momentum conser- 
vation. This leads to a problem in cases when some of the generated photon momenta are 
omitted. 

In general, it is solved by applying an unique description of how to evaluate (3 functions. 
Two approaches have been developed in the past to deal with this situation, namely the 
projection method [3] and the extrapolation method [4]. In the following we will briefly 
discuss both methods in order to point out their advantages and disadvantages. Our own 
solution will be presented afterwards. 



reads: 
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ourselves to the description of the process e + e 



ff, i.e. a process with two incoming 
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The projection method defines for each (3 n function a mapping of the four-momenta from 
the full phase space with iV 7 photons to the restricted phase space with only n photons. 
Therefore, _/V 7 — n momenta are completely discarded. The four-momenta of the remaining 
particles are transformed such that they again comply with four-momentum conservation. 
In general, this is achieved by a sequence of scalings and boosts applied in individual 
combinations on the incoming and outgoing particles. For instance, for the leading order 
term (5q the projected momenta p[,q'j are obtained by boosting the incoming momenta p\ 
and P2 into their common cm. frame and scaling them to fulfill the condition (p^ + p' 2 ) 2 = 
(Q1 + Q2) 2 = s'. The outgoing momenta q[ and q' 2 are thus defined by boosting the original 
qj into their common cm. frame. A projection for higher order corrections is more involved 
and will be discussed at the end of this section in the framework of our new approach. 
The extrapolation method, on the other hand, defines the matrix element in any phase space 
with additional photons, i.e. the momenta are kept untouched but the matrix element is 
rewritten in a suitable set of variables like e.g. angles and energies. Obviously, there is no 
universal way to construct this transformation and therefore this has to be done for each 
matrix element separately. 

As an example we briefly review the corresponding algorithm in the program KoralZ, that 
is dealing with processes of the type e + e~ — > // with / / e~. The leading order term f3o 
is written as 



Po(Pi,qj-,h, ■ ■ -,k Nj ) = \ 



— ( S ,cos^) + — ( S ,costf 2 ) 



(2.8) 



where 



and 



costf 1 = ^- = -^ (2.9) 

cos §2 = . - . , - . = - .^ (2.10) 

are the angles between equal and opposite charged particles, respectively. In the presence of 
additional photons these angles differ, but (3q is well calculable for any number of photons, 
since the variables s' = (q± + q2) 2 , cosi?i, and cos"i?2 are defined in all cases. For the 
calculation of (3\ the variable s = (p\ + P2) 2 and the polar angle of the photon cosi? 7 are 
used in addition to cosi?i, cost?2 and s'. 

However, any approach to define matrix element corrections has to comply with two con- 
ditions: 

• First, the method has to guarantee that the original result is retained in case the 
neglected "spectator" photons have total zero four-momentum. In particular, this 
must hold true in the case of zero spectators. 

• Second, the procedure should not give rise to new large logarithms at any order 
beyond the perturbative order up to which the weight is defined. 

At this stage we would like to present a new method to define the (3 functions. We aim 
at a general approach that takes full advantage of the automatic calculation of tree-level 
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matrix elements in AMEGIC++. We do not want to interfere with the automatic calculation 
as such, which relies on the construction of four-momenta and their explicit conservation. 
Hence, the method of choice is rather a realization of the projection method. 
Remember that the first order infrared finite correction term is 

Pi(Pi,Qj,k n ) = Pi~\pi,qj,k n ) - S(pi,k n ) 0o\pi,qj) ■ (2.11) 

For this correction term the projection method is well defined as long as there is exactly 
one additional photon. Then we can apply the PM as introduced for the leading order 
term (3q by redefining j3\ as 

Pi{Pi, qj, K) = Pi\pi, qj, K) - S{pi, k n ) ffl (p' i: q'j) . (2.12) 

For more than one photon the arguments of the differential cross section have to be 
projected as well. Naively, one might want to apply the same technique as for (3q with the 
photon ki included in the outgoing momenta, but it turns out that this is not sufficient to 
describe events with additional hard photons. 
Our solution to this problem is the following: 

1. The second term, consisting of the eikonal factor multiplied with the Born term cross 
section is still treated as described above. Here, we observe that 

• the scale of the evaluation of the eikonal is s, since it is evaluated with the 
original vectors p\ and p2- This means that the effect of all other photons is 
ignored which is yields the correct behavior in order to cancel the eikonal factors 
in the YFS exponent. 

• The scale of the Born cross section, i.e. the scale of the s-channel propagator, 
is s'. 

2. In the first term we have to obtain a similar structure. This is achieved by the 
following steps: 

• To calculate the photon part at the scale s, it is necessary to leave the incom- 
ing momenta as well as the photon momentum untouched, but the outgoing 
momenta have to be modified. The new momenta q\ and fa are fixed by the 
constraint 

Pi + P2 = qi + fa + k n . 
Thus, the differential cross section with one photon is given by 

P < i ) (Pi,P2,qi,fa,k n ) . (2.13) 

• In the previous step, the scale of the s-channel propagator has been changed to 
s = (qi + fa) 2 which can modify the cross section quite drastically. We cure this 
by a correction factor 

c = Po(p'i,P2,qi,Q2) ( 2 M ) 
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Thus, the right peak structure is restored. The new incoming momenta are 
naturally defined via s = (pi +P2) 2 . 

Altogether, p\ is defined as 

Po Cp' 9' ) 

Pi(Pi,qj,k n ) = v J. pi(pi,qj,kn) ~ S(pi,k n ) Po(Pi,qj) ■ (2.15) 
PO \Qii Qj ) 

3. Implementation 

In this section we want to discuss some aspects of the interplay of AMEGIC++ with its ISR 
module. A more thorough description of AMEGIC++ is available in the manual [2], here we 
just list the basic steps to calculate cross sections: 

1. As a first step, AMEGIC++ sets up the model in which the calculation is going to take 
place. This includes the definition of particles, the initialization of Feynman rules in 
terms of vertices and - possibly - the determination of the widths of unstable parti- 
cles. This last step already relies on the construction and evaluation of corresponding 
Feynman diagrams and their integration. 

2. The construction of Feynman diagrams is performed in three steps: 

• The construction of "empty" topologies, i.e. without the specification of lines 
and vertices, 

• their mapping on "filled" Feynman diagrams with specified propagators and 
vertices, 

• and, finally, their translation into helicity amplitudes which can be stored in 
library files. 

3. The integration of the diagrams is done with help of an adaptive multichannel method 
[5] . The channels are constructed by inspection of the kinematical structure of cor- 
responding Feynman diagrams. 

Implementing ISR in the framework of YFS affects the last two steps. First, AMEGIC++ 
has to provide one additional matrix element to calculate the real photon matrix element 
corrections. The second step is therfore extended as follows: 

• An "ISR-Photon" is defined, i.e. an additional particle type is introduced, with all 
properties of an ordinary photon but exclusively coupling to the charged particles in 
the initial state. 

• One "ISR-Photon" is added to the list of outgoing particles. Correspondingly, 
"empty" topologies with one additional leg are constructed. 

• The Feynman diagrams for this additional process are generated and translated into 
helicity amplitudes as pointed out above. 
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During the integration, step 3, the following changes take place: 

• The reduced cm. energy \fs' is determined, and a number of ISR photons is generated 
along the lines of the algorithm in Appendix A. With this energy a set of incoming 
momenta is defined in their cm. frame. 

• The momenta of the final state are generated and the leading order cross section is 
calculated with a multichannel method. So far the effect of ISR is only manifest by 
the reduced cm. energy. 

• The ISR weight is calculated using the generated final state momenta. Therefore, 
the matrix element prepared in step 2 with one additional photon as well as the 
Born matrix element are used according to our method introduced in section 2. The 
obtained weight is combined with the phase space weight. 

• Finally, a boost of all incoming and outgoing momenta is performed in a way that all 
momenta, i.e. including the ISR photons, comply with four-momentum conservation. 
Note that this step does not modify the result of the integrated cross section but it 
is mandatory, when AMEGIC++ is used inside a Monte Carlo generator. 

4. Results 

To demonstrate the quality of our approach we will discuss two examples: 

• Processes of the type e + e~ — > 7* , Z* — > // at LEP2 energies including comparison 
with both experimental data and other generators (KoralZ). 

• Processes of the type — > j*,Z*,h* — > // at higher energies including the 
effect of propagating Higgs bosons, both from the Standard Model and its minimal 
supersymmetric extension (MSSM). 

e+e" - // 

The effect of ISR on processes of the type e + e~ — ► // at LEP2 energies is tremendous. 
Due to the "Radiative Return" huge numbers of events actually reside on the Z° resonance 
instead of the beam energy. 

The effect is twofold: First, the visible energy in the detector tends to be reduced to 
energies of the order of the Z mass accompanied by an increase of the total cross section. 
Secondly, the forward-backward asymmetry due to the axial coupling of the Z boson tends 
to be shifted correspondingly. 

We will highlight this effect in a series of figures. Fig. 1 exhibits the total cross section 
for e + e~ — ► /U + /z _ (7) in comparison with data from the DELPHI collaboration [6]. Fig. 2 
depicts the combined forward-backward asymmetry in fj, and r pair production at LEP2 
compared to data taken by OPAL [7]. Finally, Fig. 3 shows the differential cross section 
with respect to the cosine of the angle between equally charged incoming and outgoing 
leptons compared to data from L3 [9]. To guide the eye, we have included the uncorrected 
Born level prediction with a black line in the first two figures. 
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In all three cases different cuts on s' have been applied, and in all cases the agreement of 
our generator with data is very good. Note that we did not take into account any FSR 
effects. 
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Figure 1: Energy dependence of the cross 
section a at energies above the Zq boson 
mass in comparison with DELPHI data [6]. 
Two different cuts on s' are applied, includ- 
ing and excluding a radiative return to the 
Z° pole. 
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Figure 2: Forward-backward asymmetry 
above the resonance of the Zq boson in com- 
parison with combined muon and tau pair 
data from OPAL [7]. Again, different cuts 
on the effective cm. energy y/p are plotted 
together with the Born expectation (without 
ISR). The open symbols mark preliminary 
data [8]. 



.Q 
Q. 



v> 6 

o 

o 

"D 

1 



"i 1 i 1 i 1 i 1 r 

^Js T >75GeV 

* e + e" -> |o. + |T 
■ e + e" -> x + x" 

-NJsVs > 0.85 

e + e" -> 
□ e + e" -> x + x" 



"i 1 i 1 i 1 r 



1 



J , I , I , I , I i I , L 



-1 -0.8 -0.6 -0.4 -0.2 -0 0.2 0.4 0.6 0.8 1 

cos e 

Figure 3: Differential cross section at y/s = 189 GeV in comparison with L3 data [9]. Different 
s' cuts have been applied. 



In Figs. 4 and 5 the distribution of the effective cm. energy ys' of the fi + fi pair and the 
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distribution of photon energies are shown in comparison with results from KoralZ (second 
order). We find remarkably good agreement over almost the whole energy range in the 
\fs' distribution. The energy distribution is dominated by the contribution of the hardest 
photon which is equally well described. The huge peak at 72.5 GeV corresponds to the 
radiation of one photon to reach the Z-resonance. Similarly, the edges at 48.9 GeV and 
36.25 GeV correspond to a Z return via emissions of two hard photons. In this region our 
result differs slightly from KoralZ, since the interference effects become more pronounced. 
However, since the first order matrix element is used for the correction of every explicitly 
generated photon (not only the hardest) some of the higher order matrix element effects 
are already included. Note that we have concentrated on ISR effects only and therefore 
switched off FSR and weak corrections in KoralZ. In case only the first order correction of 
KoralZ is used, the results would be nearly indistinguishable. 
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Figure 4: Distribution of the effective cm. 
energy \f~P of muon pairs at a beam energy 
of ^fs = 189 GeV. First order results by 
AMEGIC++ are compared with the second or- 
der result of KoralZ. 
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Figure 5: Energy distribution of ISR pho- 
tons in comparison with KoralZ. The pho- 
tons are ordered by their energies. Again, 
solid lines are results produced by AMEGIC++, 
dotted lines are results from KoralZ. 



Muon Collider 

As another application we consider a possible muon collider [10]. Its main purpose will 
be the precise measurement of the properties of the Higgs boson once it is found. The 
muons are used for two reasons: First, lepton colliders provide comparably pure initial 
states allowing for spectroscopical measurements (like for instance the exact determination 
of branching ratios of rare decays like h — * 77). Second, muons have a non- negligible 
coupling to the Higgs boson which in turn will be produced as an s-channel resonance. 
In Fig. 6 and Fig. 7 the total cross section and the forward-backward asymmetry for 
— * bb are plotted with and without ISR effects assuming a Standard Model Higgs 
boson with a mass of 115 GeV. The small width of roughly 4 MeV is an experimental 
challenge, e.g. effects of beam resolution have to be considered. The main effect of ISR is 
to increase the background due to radiative return events to the Z° peak. However, this 
background can be reduced via an appropriate s' cut. Of course, a radiative return to the 
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Higgs peak is also possible, but since the cross section is tiny the effect is small and only 
visible close to the peak. 
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Figure 6: Energy dependence of the cross 
section a around the resonance of the /ism 
boson in the process /z + pT — > bb. 
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Figure 7: Forward-backward asymmetry 
around the resonance of the /ism boson in 
the process fi + pT — > bb. 



Since our method is general, we can discuss further the extension to any scalar boson 
propagating in the s-channel. As an example we consider the case of Higgs bosons in the 
Minimal Supersymmetric extension of the Standard Model (MSSM) [11, 12]. In this model 
there are two more neutral Higgs bosons which tend to be considerably heavier than the 
light one. For the following discussion let us constrain ourselves to a supergravity inspired 
supersymmetry breaking scenario [13, 16] with parameters to be found in the Appendix 
B. In this scenario we have a mass of m# = 395.0 GeV for the CP-even and and a mass 
of uiaq = 394.5 GeV for the CP-odd Higgs boson. The total cross section can be found 
in Fig. 8 and the influence on the forward-backward asymmetry is given in Fig. 9. At this 
point in the mSUGRA parameter space the individual peaks can not be resolved and the 
two Higgs bosons appear as one single, comparably broad, resonance only. The effect of 
ISR is a slight shift of the peak to larger cm. energies and an increase of background due 
to radiative return to the lower lying Z boson resonance. However, the latter effect can 
again be reduced by a suitable cut on the effective cm. energy s' . 

In striking contrast to the case of the lighter SM Higgs boson, in this case the radiative 
return to the peak structure from higher energies is important, see Figs. 10 and 11. Note 
that the contribution of a photon, a Z° boson, and all three Higgs bosons propagating in the 
s-channel as well as all interference terms have been taken into account. The effect of the 
interplay of the different terms becomes apparent in angular distributions, cf. Fig. 12, where 
different cuts on the cm. energy have been applied. The Born term is dominated by the 
off-peak Z boson contribution and exhibits therefore a slight angular dependence. When 
taking the radiative return to the heavy Higgs bosons peak into account, the distribution 
is only shifted due to the angular blind scalar coupling structure. The predominating 
contribution of the return to the Z resonance restores the angular distribution already 
known from Fig. 3. 
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Figure 8: Energy dependence of the cross 
section a around the resonance of the heavy 
scalar Higgs boson Hq and the pseudo scalar 
Higgs boson A in the process — * bb 
for different s' cuts, as well as without ISR. 
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Figure 9: Asymmetry around the resonance 
of the heavy scalar Higgs boson Hq and the 
pseudo scalar Higgs boson Aq with and with- 
out ISR. 




Figure 10: Distribution of the effective cm. Figure 11: Detail of Fig. 10. 

energy vs 7 for the process — > bb at a 

beam energy of £?cms=400 GeV. 



5. Conclusion and Outlook 

In this paper we presented an implementation of the successful YFS algorithm into the 
matrix element generator AMEGIC++. To make full use of the possibilities offered by au- 
tomatic matrix element generation with help of the helicity method and fully automatic 
integration of the matrix element we developed a new method to deal with the first order 
matrix element correction to multi photon emissions. We proved the quality of our method 
by two examples: 

1. For lepton pair production at energies at and above the Z boson resonance, we 
have confronted our simulator with experimental data from the LEP experiment and 
another ISR generator, KoralZ. In both cases we found encouraging agreement. 
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Figure 12: Angular distribution of the b-quark with respect to the incoming muon. Different 
cuts applied on y/s' in order to include return events to the Z° or to the heavy Higgs bosons Hq 
and Aq only. For comparison the Born expectation is plotted too. 

2. The application of our approach to the s-channel production of Higgs bosons at 
a muon collider has illustrated the generality of our approach. Let us note that 
extensions to other bosons in the s-channel with other couplings is straightforward 
and will be addressed in the near future. 

3. Furthermore, the interface of AMEGIC++ with the parton shower implemented in 
APACIC++ [14]provides a powerful tool for the simulation of hadronic final states 
including the effect of both ISR and the final state parton shower. 

We feel obliged to comment briefly on possible improvements: 

• Definitely, effects of Final State Radiation and its interference with ISR are missing 
and should be implemented. Additionally, some types of one-loop corrections like 
QCD corrections in the final state, electroweak corrections or box diagrams have not 
been implemented yet . Also, an extension to second order corrections should be 
done to achieve the precision needed for a proper interpretation of the anticipated 
precision data. 

• Apart from s-channel processes, at e + e~ colliders processes like ZZ production play 
a crucial role and deserve a similar treatment. In the case of for instance WW 
production an extension of our method is not a trivial task. 

1 During the time of writing we became aware of [15], where an ISR calculation for a Higgs boson 
production at a muon collider was considered. This was based on the structure function approach and 
already includes FSR, IF-interference and QCD effects. Their findings are consistent with our Monte Carlo 
results. 



Born 

Nl? > 390 GeV 
\ls" > 40 GeV 
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A. The MC algorithm 

In this appendix we would like to discuss the Monte Carlo procedure for the generation of 
ISR photons, further details of the first three steps can be found in [3]. 

1. Choice of the reduced cm. energy squared 

First, the reduced cm. energy squared is determined. This is done according to a 
"crude distribution" that consists of a fast approximation for the energy dependence 
of the cross section. This approximation already exhibits the correct peak structure as 
given by the s-channel propagators. Additionally, it includes the leading logarithmic 
photon factor. The corresponding probability density in terms of the energy fraction 
v = 1 — s'/s reads 

f(v) = \^r— V(v) (A.l) 

with 

9ry is \ 

= (l-i/)<7Bom(s') J{v) and 7 = — ( In - 1 In e . (A.2) 

The energy fraction v can be determined according to f{y) by standard Monte Carlo 
sampling techniques. In the equation above J(y) is an approximate form of the 
Jacobian that will be discussed in step 3(b). 

2. Choice of number of resolvable photons 

Having determined the effective center of mass energy squared s' = s(l — v), we 
proceed by choosing the number of resolvable photons which in turn have to be 
generated explicitly. In the logarithmic approximation the photon number iV 7 follows 
a shifted Poisson distribution: 

w " l) - ^y^ry. ■ <A ' 3) 

Obviously, its mean value depends on both s' and the photon resolution parameter 
e. It is given by 

= ??L]n(— y \]n(-)Q(v-e). (A.5) 
7T \m e J \£J 
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This can be identified in the master formula, Eq.(2.5), considering the additional 
constraint on the photon energies (the second 0-function) not to exceed s'. 

3. Construction of four-momenta 

The construction of four-momenta for the photons is achieved in two steps: 

(a) Angles: The photon angles are now individually appointed to each photon with 
the help of the eikonal S. Clearly, S exhibits no dependence on the azimuthal an- 
gle <p, hence this angle is distributed uniformly. In contrast, S depends strongly 
on the polar angle which in turn is chosen according to 



I / 4m 2 

/(COS t?) = 7 nw n ST With P = \ 1 • (A.6) 

y 1 (1 - COS +0 COST?) ' V S V ' 

Note that terms proportional to the lepton mass have been discarded at this 
stage. Therefore, a correction weight is introduced for each photon 

™an g leW = , (A.7) 

where the exact distribution /(#) is given by 

m 2 ( 1 1 \ 

m = /(cos *) - [ {1 _ pcos ^ + (1 + /?cos ^ )2 J • (A.8) 

Thus, the weight to reproduce the correct mass-dependence reads 

W mass = JJ^angle(^i) • (A.9) 

i=l 

(b) Energies: The next step is the assignment of energies to each of the photons. 
In case that exactly one photon above the threshold has been generated, the 
energy of this photon is already completely constrained by v and the chosen 
photon angles due to overall four-momentum conservation. In the case that 
more than one photon have been generated, the additional photon energies are 
chosen according to dk°/k°, i.e. the eikonal S. Then a scaling of all photons 
is performed in order to achieve four-momentum conservation. Thereby an 
Jacobian 

w ' ) -K 1+ 7ra) <A10) 

with the abbreviations 

K 2 P 2 N ~* 

A = WW 1 K = ^ kn (A - n) 

is introduced. It had been taken into account already in step 1 as J assuming 
A = 1, since the photon momenta have not been know at this time. A correction 
weight 

^jac = ~ . , (A. 12) 
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is introduced to cure this simplification. Furthermore, all photon energies have 
to be above the limit after rescaling, in order to avoid double counting. 

This leads to an additional weight 




4. Matrix element correction 

The next step is the calculation of the matrix element correction. This is done by 
introducing an additional weight according to 



5. Re- weighting 

The last step is the calculation of the final weight in order to cure the approximations 
made during the choice of the reduced cm. energy (step 1) and the determination of 
the photon momenta (step 3) as well as to incorporate the matrix element corrections 



This weight can now be used in a rejection method to produce unweighted events. 

The ISR simulation in AMEGIC++ differs from this general description mainly in steps 1 and 
4. 

In step 1, a multichannel method is used to determine s' to incorporate complicated peak 
structures. It can therefore be understood as the first part of the multichannel determina- 
tion of the final state momenta. A library of fast // — > /'/' cross sections for arbitrary 
scalar and vector s-channel resonances enables the simulation of many processes in the 
SM and beyond. Due to the rejection mechanism these cross sections can be also used for 
multi-jet productions. 

In step 4 the evaluation of the matrix element correction is modified as described in section 
2, thereby taking full advantage of the features of AMEGIC++. 

B. Program Parameters 

The parameters used for the calculation of results presented in section 4 are given in Tab. 1 



To include higher order effects in the Yukawa coupling we have used a running b-quark 
mass according to the well known relation 




(A.14) 



WMG = W m3ss W ]3JC W £ WME ■ 



(A.15) 



and Tab. 2. 




(B.l) 



with 



7m,0 = 1 , A) 




(B.2) 



«QED 


1./137.035995 


sin 2 d w 


0.23117 


M z 


91.1882 GeV 


Tz 


2.4952 GeV 


m b 


5.0 GeV 


m c 


1.3 GeV 




105.658 MeV 


m T 


1.777 GeV 


V 


246.2 GeV 


m h am 


115.0 GeV 




3.8 MeV 



Table 1: Standard Model parameter used 
for the simulation. 



M 


100. 


GeV 


m Ho 


395.0 GeV 


M V 2 


250. 


GeV 


r H 


1.07 GeV 


A 


-100. 


GeV 


m Ao 


394.5 GeV 


tan/3 


10. 






1.29 GeV 


A* 


+1 








m t 


175. 


GeV 







Table 2: The mSugra point SPS1 according 
to [16] and the corresponding properties of 
the heavy Higgs bosons in the MSSM. The 
widths of the Higgs bosons have been calcu- 
lated with AMEGIC++. 



- 18 - 



References 



D. R. Yennie, S. C. Frautschi, and H. Suura, Ann. Phys. 13 (1961), 379-452. 

F. Krauss, R. Kuhn and G. Soff, JHEP 0202 (2002) 044 

S. Jadach and B. F. L. Ward, Comput. Phys. Commun. 56 (1990), 351-384. 

S. Jadach, B. F. Ward and Z. Was, Comput. Phys. Commun. 79 (1994) 503. 

R. Klciss and R. Pittau, Comput. Phys. Commun. 83 (1994) 141. F. A. Berends, R. Pittau 
and R. Klciss, Nucl. Phys. B 424 (1994) 308. 

P. Abreu et al. [DELPHI Collaboration], Phys. Lett. B 485 (2000) 45. P. Abrcu et al. 
[DELPHI Collaboration], Eur. Phys. J. C 11 (1999) 383. P. Abreu et al. [DELPHI 
Collaboration], Eur. Phys. J. C 16 (2000) 371. 

G. Abbicndi et al. [OPAL Collaboration], Eur. Phys. J. C 13 (2000) 553. 

OPAL Collab., OPAL Physics Note PN469 (February 2001) and relerences therein. 

M. Acciarri et al. [L3 Collaboration], Phys. Lett. B 479 (2000) 101. 

C. M. Ankcnbrandt et al., Phys. Rev. ST Accel. Beams 2 (1999) 081001. 

see for instance: H. E. Haber and G. L. Kane, Phys. Rept. 117 (1985) 75. 

V. D. Barger, M. S. Bergcr, J. F. Gunion and T. Han, Phys. Rept. 286 (1997) 1. 

see for instance: H. P. Nillcs, Phys. Rept. 110 (1984) 1. 

R. Kuhn, F. Krauss, B. Ivanyi and G. Soff, Comput. Phys. Commun. 134 (2001) 223. 
S. Dittmaier and A. Kaiser, hep-ph/0203120. 

B. C. Allanach et al, in Proc. of the APS/DPF/DPB Summer Study on the Future of 
Particle Physics (Snowmass 2001) ed. R. Davidson and C. Quigg. 



- 19 - 



